\documentclass[a4paper,12pt]{article}
\usepackage{graphicx}
\usepackage{listings}
\usepackage[justification=centering]{caption}
\usepackage{subfigure}  
\usepackage{multirow}
\usepackage{balance}
\usepackage{gensymb}
\lstset{language=Matlab}
\lstset{breaklines}
\lstset{extendedchars=false}

\usepackage{amsmath,amsfonts,amsthm,amssymb}
\theoremstyle{definition}
\newtheorem{thm}{Theorem}
\newtheorem{exmp}{Example}
\newtheorem{defn}{Definition}
\newtheorem{lema}{Lemma}
\newtheorem{prop}{Proposition}
\newtheorem{coro}{Corollary}

\newcommand{\GreekAlpha}[1]{\uppercase\expandafter{\romannumeral#1}}
\newcommand{\Greekalpha}[1]{\lowercase\expandafter{\romannumeral#1}}
\renewcommand{\baselinestretch}{1.25}

%------------setlength----------------%
\setlength{\textwidth}{162mm}
%\setlength{\textheight}{190mm}
\setlength{\textheight}{231mm}
\setlength{\headheight}{-0.1cm}
\setlength{\topmargin}{-0.1cm}
\setlength{\oddsidemargin}{-0cm}
\setlength{\evensidemargin}{-0cm}

\setlength\columnsep{12pt}
\setlength\columnseprule{0.4pt}

\setlength{\parskip}{1mm}
\setlength{\unitlength}{1mm}
\setlength{\parindent}{2em}

\title{Numerical Analysis Homework \#6}

\author{Li Zhiqi\quad 3180103041 , }

\begin{document}
\maketitle
\section{Theoretical questions}
\GreekAlpha{1}.\\
\indent (a) We consider the interpolation polynomial satisfying the condition
required. The table of divided differences is as follows.
\begin{table}[!htp]
	\centering
	\setlength{\tabcolsep}{0.6mm}{
	\begin{tabular}{c|c c c c}
		$-1\ $ & $f(-1)\ $ &  &  &\\	
		$0\ $ & $f(0)\ $ & $f(0)-f(-1)\ $ & &  \\			
                $0\ $ & $f(0)\ $ & $f'(0)\ $ & $f'(0)-f(0)+f(-1)\ $ &  \\
                $1\ $ & $f(1)\ $& $f(1)-f(0)\ $ & $f(1)-f(0)-f'(0)\ $ & $\frac{f(1)-2f'(0)-f(-1)}{2}\ $
	\end{tabular}
}
\end{table}\\
Hence
  \begin{align*}
    &p_3(y;-1,0,0,1;t) \\
    &= f(-1)+(f(0)-f(-1))(t+1)+(f'(0)-f(0)+f(-1))t(t+1) \\
    &+ \frac{f(1)-2f'(0)-f(-1)}{2}t^2(t+1).\\
  \end{align*}
  Then we have
  \begin{align*}
    &\int_{-1}^1p_3(y;-1,0,0,1;t)\text{d}t\\
    &=\int_{-1}^1(f(-1)+(f(0)-f(-1))+(f'(0)-f(0)+f(-1))t^2+\frac{f(1)-2f'(0)-f(-1)}{2}t^2)\text{d}t\\
    &=2f(0)+\frac{2}{3}[f'(0)-f(0)+f(-1)+\frac{f(1)-2f'(0)-f(-1)}{2}]\\
    &=\frac{1}{3}[f(-1)+4f(0)+f(1)] = \int_{-1}^1y(t)\text{d}t - E^S(y),
  \end{align*}
  which completes the proof.\\
  
  (b) By (a) and Theorem 3.33,
  \begin{align*}
    &E^S(y) = \int_{-1}^1(y(t) - p_3(y;-1,0,0,1;t))\text{d}t\\
    &=\int_{-1}^1\frac{f^{(4)}(\xi)}{4!}(t+1)t^2(t-1)\text{d}t\\
    &=\frac{f^{(4)}(\xi)}{24}\int_{-1}^1(t^4-t^2)\text{d}t = -\frac{f^{(4)}(\xi)}{90},
  \end{align*}
  where $\xi \in (-1,1)$.\\

  (c) Firstly, we change the integration interval. Set
  $x = a + (b-a)t$, then
  \begin{align*}
    \int_b^ay(x)\text{d}x &= \int_{-1}^1y(a+(b-a)t)(b-a)\text{d}t\\
                          &=(b-a)\int_{-1}^1p_3(y';-1,0,0,1;t)\text{d}t + E^S\\
    &=\int_a^bp_3(y;a,\frac{a+b}{2},\frac{a+b}{2},b;x)\text{d}x +
      E^s\\
    &=\frac{b-a}{6}[f(a)+4f(\frac{a+b}{2})+f(b)] + E^s.
  \end{align*}
  Hence
  \begin{align*}
    E^S &=
    \int_a^b(y(x)-p_3(y;a,\frac{a+b}{2},\frac{a+b}{2},b;x))\text{d}x\\
        &=\int_{a}^b\frac{f^{(4)}(\xi)}{4!}(x-a)(x-\frac{a+b}{2})^2(x-b)\text{d}x\\
    &=-\frac{(b-a)^5}{2880}f^{(4)}(\xi),
  \end{align*}
  where $\xi \in (a,b)$.\\
  Next we use above conclusion in composite Simpson's rule
  \begin{align*}
    I_n^S(f) =
    \frac{h}{3}[f(x_0)+4f(x_1)+2f(x_2)+4f(x_3)+2f(x_4)+\cdots+4f(x_{n-1})+f(x_n)]
  \end{align*}
  by subintervals, and sum up the errors, we have
  \begin{align*}
    E_n^s(f) &= -\sum_{i=1}^{n/2}\frac{(2h)^5}{2880}f^{(4)}(\xi_i)\\
    &= -\frac{b-a}{180}h^4[\frac{2}{n}\sum_{i=1}^{n/2}f^{(4)}(\xi_i)]
      = -\frac{b-a}{180}h^4f^{(4)}(\xi),
  \end{align*}
  where the last step comes from the intermediate value Theorem C.39
  and the fact $f^{(4)} \in \mathcal{C}(a,b)$. Here $\xi \in (a,b)$.\\
  
  
  \noindent \GreekAlpha{2}.\\
  \indent (a) By Theorem 6.18, we need
  \begin{align*}
    |E_n^T(f)| = |\frac{b-a}{12}h^2f''(\xi)| =
    \frac{1}{12n^2}|f''(\xi)| \le 0.5 \times 10^{-6}.
  \end{align*}
  Notice $f''(x) = (4x^2-2)e^{-x^2} < 2e^{-1}, \forall x\in (0,1)$, hence
  \begin{align*}
    n \ge \sqrt{\frac{2e^{-1}}{12\times0.5 \times 10^{-6}}}
    \approx 350.18.
  \end{align*}
  Then 351 subintervals are required.\\

  (b) Similarily, by Theorem 6.20, we need
  \begin{align*}
    |E_n^S(f)| = |\frac{b-a}{180}h^4f^{4}(\xi)| =
    \frac{1}{180n^4}|f^{(4)}(\xi)| \le 0.5 \times 10^{-6}.
  \end{align*}
  Notice $f^{(4)}(x) = (16x^4-48x^2+12)e^{-x^2} < 12, \forall x\in
  (0,1)$, hence
  \begin{align*}
    n \ge (\frac{12}{180\times0.5 \times 10^{-6}})^{\frac{1}{4}}
    \approx 19.11.
  \end{align*}
  Then 20 subintervals are required.\\

  \noindent \GreekAlpha{3}.\\
  (a) We have
  \begin{align*}
    \left \langle 1,\pi_2 \right \rangle &= \int_0^{+\infty}1\cdot
    \pi_2(t)\rho(t)\text{d}t =
    \int_0^{+\infty}(t^2+ at + b)e^{-t}\text{d}t = 2+a+b = 0;\\
    \left \langle t,\pi_2 \right \rangle &= \int_0^{+\infty}t\cdot
    \pi_2(t)\rho(t)\text{d}t =
    \int_0^{+\infty}t(t^2+ at + b)e^{-t}\text{d}t = 6+2a+b = 0,
  \end{align*}
  which implies $a = -4, b = 2$. Hence
  \begin{align*}
    \pi_2(t) = t^2-4t+2.
  \end{align*}

  (b) The roots of $\pi_2(t)$ are $t_1 = 2-\sqrt{2},t_2 = 2+\sqrt{2}$.
  To get $\omega_1$ and $\omega_2$, consider
  \begin{align*}
    \int_0^{+\infty}e^{-t}\text{d}t &= 1 = \omega_1 + \omega_2,\\
    \int_0^{+\infty}te^{-t}\text{d}t &= 1 = t_1\omega_1 + t_2\omega_2,
  \end{align*}
  which implies $\omega_1 = \frac{2+\sqrt{2}}{4},\omega_2 =
  \frac{2-\sqrt{2}}{4}$. Hence we have
  \begin{align*}
    \int_0^{+\infty}f(t)e^{-t}\text{d}t = \frac{2+\sqrt{2}}{4}f(2-\sqrt{2}) +
    \frac{2-\sqrt{2}}{4}f(2+\sqrt{2}) + E_2(f),
  \end{align*}
  with
  \begin{align*}
    E_2(f) =
    \frac{f^{(4)}(\tau)}{4!}\int_0^{+\infty}e^{-t}(t^2-4t+2)^2\text{d}t
    = \frac{f^{(4)}(\tau)}{6},
  \end{align*}
  where $\tau \in (0,\infty)$, and the first equation follows from
  Theorem 6.36.

  (c) With $f(t) = \frac{1}{1+t}$, we compute
  \begin{align*}
    \int_0^{+\infty}f(t)e^{-t}\text{d}t = \frac{2+\sqrt{2}}{4}f(2-\sqrt{2}) +
    \frac{2-\sqrt{2}}{4}f(2+\sqrt{2}) \approx 0.571428571, 
  \end{align*}
  combining with $f^{(4)}(t) = \frac{24}{(1+t)^5} < 24, \forall t \in (0,\infty) $,
  \begin{align*}
    E_2(f) = \frac{f^{(4)}(\tau)}{6} < 4.
  \end{align*}
  By using the exact value $I = 0.596347361\cdots$, we have
  \begin{align*}
    E_2^t(f) = I - \frac{2+\sqrt{2}}{4}f(2-\sqrt{2}) +
    \frac{2-\sqrt{2}}{4}f(2+\sqrt{2}) \approx 0.024918790 \ll 4,
  \end{align*}
  and
  \begin{align*}
    E_2(f) = \frac{f^{(4)}(\tau)}{6} = \frac{4}{(1+\tau)^5} = E_2^t(f)\approx 0.024918790,
  \end{align*}
  which implies $\tau \approx 1.7612556$.\\

  \noindent \GreekAlpha{4}.\\
  (a) Notice that
  \begin{align*}
    p(x_i) &= \sum_{m=1}^n[h_m(x_i)f_m+q_m(x_i)f'_m] = (a_i+b_ix_i)f_i
    + (c_i+d_ix_i)f'_i = f_i,\\
    p'(x_i)&=
             \sum_{m=1}^n[h'_m(x_i)f_m+q'_m(x_i)f'_m]=[b_i+(a_i+b_ix_i)2l'_i(x_i)]f_i+[d_i+(c_i+d_ix_i)2l'_i(x_i)]f'_i
             = f'_i.
  \end{align*}
  Set
  \begin{align*}
    \left\{\begin{matrix}
  &a_i+b_ix_i=  1 \\
  &c_i+d_ix_i=  0 \\
  &b_i+(a_i+b_ix_i)2l'_i(x_i)=  0 \\
  &d_i+(c_i+d_ix_i)2l'_i(x_i)=  1
\end{matrix}\right. \qquad \forall i = 1,2,\cdots,n.
  \end{align*}
  Then we get
  \begin{align*}
    a_m = 1+2x_ml'_m(x_m), b_m = -2l'_m(x_m), c_m = -x_m, d_m = 1, \forall m = 1,2,\cdots,n.
  \end{align*}

  (b) By (a), $\forall f \in \mathbb{P}_{2n-1}$, $p(f) = f$. Hence we consider
  \begin{align*}
    I_n(f) = I_n(p(f)) &= \int \rho\sum_{k=1}^n(h_kf_k+q_kf'_k)\\
              &=\sum_{k=1}^n[f(x_k)\int \rho h_k + f'(x_k)\int \rho q_k]\\
    &=\sum_{k=1}^n[\omega_kf(x_k) + \mu_kf'(x_k)].
  \end{align*}
  Here
  \begin{align*}
    \omega_k &= \int \rho(t)h_k(t)\text{d}t = \int \rho(t)(1+2x_kl'_k(x_k)-2l'_k(x_k)t)l_k^2(t)\text{d}t, \\
    \mu_k &= \int \rho(t)q_k(t)\text{d}t = \int \rho(t)(t-x_k)l_k^2(t)\text{d}t.
  \end{align*}

  (c) Notice
  \begin{align*}
    \mu_k &= \int \rho(t)(t-x_k)l_k^2(t)\text{d}t = 0\\
    &\Rightarrow \int \rho(t)v_n(t)l_k(t)\text{d}t = 0,\quad \forall k = 1,2,\cdots,n.
  \end{align*}
  Combining the fact that $\text{span}\{l_1(t),l_2(t),\cdots,l_n(t)\}
  = \text{span}\{1,t,\cdots,t^{n-1}\}$, we get
  \begin{align*}
    \int \rho(t)v_n(t)p(t)\text{d}t = 0,\quad \forall p \in \mathbb{P}_{n-1}.
  \end{align*}

  \noindent \GreekAlpha{5}.\\
  We denote the actual input value with random errors $\epsilon$ by
  $\tilde{u}$, and the exact value by $u$. Then
  \begin{align*}
    |u''(\bar{x}) - D^2u(\bar{x})| & \le |u''(\bar{x}) -
                                     D_T^2u(\bar{x})| +
                                     |D_T^2u(\bar{x}) -
                                     D^2u(\bar{x})|.
  \end{align*}
  The first term can be controled by Taylor formula:
  \begin{eqnarray}
    \footnotesize
  \begin{aligned}
    &|u''(\bar{x}) - D_T^2u(\bar{x})| = |u''(\bar{x}) -
                                       \frac{u(\bar{x}-h)-2u(\bar{x})+u(\bar{x}+h)}{h^2}|\\
    &\!=\!|u''(\bar{x})\!
    -\!\frac{u(\bar{x})\!-\!hu'(\bar{x})\!+\!\frac{h^2}{2}u''(\bar{x})\!-\!\frac{h^3}{6}u'''(\bar{x})\!+\!\frac{h^4}{24}u^{(4)}(\xi_1)\!-\!2u(\bar{x})\!+\!u(\bar{x})\!+\!hu'(\bar{x})\!+\!\frac{h^2}{2}u''(\bar{x})\!+\!\frac{h^3}{6}u'''(\bar{x})\!+\!\frac{h^4}{24}u^{(4)}(\xi_2)}{h^2}|
    \\
    &=\frac{h^2}{12}|u^{(4)}(\xi)|, \nonumber
  \end{aligned}
  \end{eqnarray}
  where $\xi \in [\bar{x}-h,\bar{x}+h]$. On the other hand, the second
  terms
  \begin{align*}
    |D_T^2u(\bar{x}) - D^2u(\bar{x})| &=
                                        |\frac{u(\bar{x}-h)-2u(\bar{x})+u(\bar{x}+h)}{h^2}
                                        -
                                        \frac{\tilde{u}(\bar{x}-h)-2\tilde{u}(\bar{x})+\tilde{u}(\bar{x}+h)}{h^2}|\\
    &\le
      \frac{1}{h^2}(|u(\bar{x}-h)-\tilde{u}(\bar{x}-h)|+2|u(\bar{x})-\tilde{u}(\bar{x})|+|u(\bar{x}+h)-\tilde{u}(\bar{x}+h)|)\\
    &\le \frac{4E}{h^2}.
  \end{align*}
  So far we have completed the proof that
  \begin{align*}
    |u''(\bar{x}) - D^2u(\bar{x})| \le \frac{h^2}{12}|u^{(4)}(\xi)| + \frac{4E}{h^2},
  \end{align*}
  where $\xi \in [\bar{x}-h,\bar{x}+h]$.\\
  Notice that
  \begin{align*}
    \frac{h^2}{12}|u^{(4)}(\xi)| + \frac{4E}{h^2} \ge
    2\sqrt{\frac{h^2}{12}|u^{(4)}(\xi)|\cdot \frac{4E}{h^2}} = 2\sqrt{\frac{E|u^{(4)}(\xi)|}{3}}.
  \end{align*}
  And $\frac{h^2}{12}|u^{(4)}(\xi)|= \frac{4E}{h^2}$ implies $h =
  (\frac{48E}{|u^{(4)}(\xi)|})^{\frac{1}{4}}$.

  Next, we design a fourth-order accurate formula based on a symmetric
  stencil as follows:
  \begin{align*}
    D^4u(\bar{x}) = \frac{-u(\bar{x}-2h)+16u(\bar{x}-h)-30u(\bar{x})+16u(\bar{x}+h)-u(\bar{x}+2h)}{12h^2}.
  \end{align*}
  By similar analysis, we get
  \begin{align*}
    |u''(\bar{x}) - D^4u(\bar{x})| \le \frac{h^4}{90}|u^{(6)}(\xi)| + \frac{16E}{3h^2},
  \end{align*}
  where $\xi \in [\bar{x}-2h,\bar{x}+2h]$, and
  \begin{align*}
    \frac{h^4}{90}|u^{(6)}(\xi)| + \frac{16E}{3h^2} &=
    \frac{h^4}{90}|u^{(6)}(\xi)| + \frac{8E}{3h^2} + \frac{8E}{3h^2}\\
    &\ge
    3\sqrt[3]{\frac{h^4}{90}|u^{(6)}(\xi)|\cdot\frac{8E}{3h^2}\cdot\frac{8E}{3h^2}}\\
    &=  4\sqrt[3]{\frac{E^2|u^{(6)}(\xi)|}{30}}.
  \end{align*}
  Here $\frac{h^4}{90}|u^{(6)}(\xi)| = \frac{8E}{3h^2}$ implies $h =
  (\frac{240E}{|u^{(6)}(\xi)|})^{\frac{1}{6}}$.\\
  We find that the error bound does not depend on $h$, which means we
  can't get specify precision by decreasing $h$. Moreover, the error
  increases as $h$ decreases while $h$ is too small. And when $h$ is
  small, the fourth-order method has no advantage over the
  second-order method.\\

\noindent \GreekAlpha{5}.(\textbf{Exercise 6.41.})\\
For $D_{+} u(\bar{x})$,
\begin{align*}
\frac{u(\bar{x}+h)-u(\bar{x})}{h}-u'(\bar{x})=
\frac{u(\bar{x})+hu'(\bar{x})+\frac{1}{2}h^2u''(\bar{x})+O(h^3)-u(\bar{x})}{h}-u'(\bar{x})
=  \frac{1}{2}hu''(\bar{x})+O(h^2),
\end{align*}
hence $D_{+} u(\bar{x})$ is first-order accurate.\\
For $D_{-} u(\bar{x})$,
\begin{align*}
\frac{u(\bar{x})-u(\bar{x}-h)}{h}-u'(\bar{x})=\frac{u(\bar{x})-(u(\bar{x})-hu'(\bar{x})+\frac{1}{2}h^2u''(\bar{x})+O(h^3))}{h}-u'(\bar{x})=  -\frac{1}{2}hu''(\bar{x})+O(h^2),
\end{align*}
hence $D_{-} u(\bar{x})$ is first-order accurate.\\
For $D_{0} u(\bar{x})$,
\begin{eqnarray}
  \small
\begin{aligned}
&\frac{u(\bar{x}+h)-u(\bar{x}-h)}{2 h}-u'(\bar{x})=\\
&\frac{u(\bar{x})+hu'(\bar{x})+\frac{1}{2}h^2u''(\bar{x})+\frac{1}{6}h^3u'''(\bar{x})+O(h^4)}{2h}-\frac{u(\bar{x})-hu'(\bar{x})+\frac{1}{2}h^2u''(\bar{x})-\frac{1}{6}h^3u'''(\bar{x})+O(h^4)}{2h}-u'(\bar{x})\\
&=  \frac{1}{6}h^2u''(\bar{x})+O(h^3),\nonumber
\end{aligned}
\end{eqnarray}
hence $D_{0} u(\bar{x})$ is second-order accurate.\\

\noindent \GreekAlpha{6}.(\textbf{Exercise 6.42.})\\
\begin{table}[!htp]
	\centering
	\setlength{\tabcolsep}{0.6mm}{
	\begin{tabular}{c|c c c}
	
		$\bar{x}$ & $u(\bar{x})$ &  &  \\	
		$\bar{x}-2h$ & $u(\bar{x}-2h)$ & $u[\bar{x},\bar{x}-2h]$ &   \\			
		$\bar{x}-h$ & $u(\bar{x}-h)$ & $u[\bar{x}-2h,\bar{x}-h]$ & $u[\bar{x},\bar{x}-2h,\bar{x}-h]$   \\	
	
	
	\end{tabular}
}
\end{table}\\
We get the quadratic polynomial:
\begin{align*}
u(x)=u(\bar{x})+u[\bar{x},\bar{x}-2h](x-\bar{x})+u[\bar{x},\bar{x}-2h,\bar{x}-h](x-\bar{x})(x-\bar{x}+2h).
\end{align*}
Then
\begin{align*}
u'(\bar{x})&=u[\bar{x},\bar{x}-2h]+2hu[\bar{x},\bar{x}-2h,\bar{x}-h]\\
&=\frac{u(\bar{x})-u(\bar{x}-2h)}{2h}+2h\cdot \frac{u(\bar{x})-2u(\bar{x}-h)+u(\bar{x}-2h)}{h^2}\\
&=\frac{3u(\bar{x})-4u(\bar{x}-h)+u(\bar{x}-2h)}{2h},
\end{align*}
which is the same as formula (6.37).
  
\end{document}
